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Abstract 


A comprehensive analysis, of an axisymmetric upsetting process has been done to study 
the effects of various process variables on process parameters and deformation The 
process variables considered are the percent reduction in height (%r), the height to 
diameter ratio (H/D) and the material properties. The process parameters studied 
are the upsetting load, the distribution of contact normal stress at t he dio-wotkpioee 
interface and the contours of effective stress. The deformation is represented by the 
deformed mesh showing the bulge profile and the contours of effective strain. Effects 
of strain hardening on the process parameters and on the deformation are studied by 
considering two materials (aluminium and steel). 

The inertial forces have been neglected. The material is assumed to be elastic- 
plastic strain hardening yielding according to von Mises criterion. The dies are assumed 
to be rigid and rough. So the process is considered under the condition of full sticking at 
the die-workpiece interface. The updated Lagrangian formulation which is convenient 
for handling geometric and material nonlinearities is used. The modified Newton- 
Raphson iterative scheme is used for solving the incremental equations. Validation 
of the model is done by comparing the results with available experimental and FEM 
results. 

It is shown that, when a full sticking condition is used, the normal contact stress 
at an edge of the interface acquires an un-realistic value of infinity. Further the phe- 
nomenon of the fold-over of the free surface is observed for %r > 20 or H/D < 1. To 
avoid these things, the interfacial slip should be incorporated in the formulation. It is 
found that the H/D ratio does not have much influence on the process parameters. It is 
also observed that strain hardening doesn’t affect the bulge profile or the distribution 
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and the values of the effective strain. 
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Nomenclature 


{a} flow vector. 

A slope of the hardening curve. 

[B l ] linear strain displacement matrix. 

[B^] non-linear strain displacement matrix. 

Cijki fourth order elasticity tensor. 

C Ep t fourth order elastic-plastic tensor. 

[C] elasticity matrix. 

[C EP ] Elastic-plastic matrix. 

D diameter of the cylindrical specimen. 
e X] Green-Lagrange strain tensor. 

E Young’s modulus. 

{/} global internal force vector. 

{/} e internal force vector for the element. 

F plastic potential. 

{F} global external force vector. 

{F} e external force vector for the element. 

H height of the cylindrical specimen. 

K hardening parameter. 

[K] Stiffness matrix. 

[Kl\ linear part of the stiffness matrix. 

[Kn l] non-linear part of the stiffness matrix. 
n hardening exponent. 

N x shape function of the element. 

p hardening variable. 

{P} load vector. 

[Q] matrix containing the shape functions of the element, 
r position vector. 

R virtual work of external forces in FE formulation. 

S t] 2nd Piola-Kirchoff stress tensor 

t time 
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u T displacement in r-direction. 

u z displacement in z-direction. 

{u} displacement vector. 

{u} e elemental displacement vector. 

V volume of the body. 

6 variation. 

A increment in a quantity. 

natural strain tensor, 
ef, elastic natural strain tensor, 

efj plastic natural strain tensor. 

e eq equivalent strain. 

e p q equivalent plastic strain. 

Arj tJ incremental non-linear strain tensor. 

A scalar multiplier in flow rules. 
v Poisson’s ratio. 

p density, 

cry Cauchy stress tensor. 

deviatoric part of Cauchy stress tensor, 
oy, Jaumann stress rate. 

a eq effective ( or equivalent ) stress. 
a 0 initial yield stress. 

o y current yield stress. 

Q t j spin tensor. 


Chapter 1 


Introduction 


In forging, the material is squeezed between two or more dies in such a way that the 
metal or alloy is plastically deformed to the desired shape and size. It is also used to 
improve the metallurgical and mechanical properties of the material. This is the oldest 
of the metal working processes known to mankind since the copper age. 

Generally the forging is performed by heating the material and then applying 
a compressive force to obtain the final shape. This is called the hot forging. The 
cold forging operation is carried out at room temperature. The forging operation is 
carried out in two different ways : drawing out and upsetting. In the drawing out 
operation, the metal gets elongated with a reduction in the cross-sectional area. For 
this purpose, the force is applied in a direction, perpendicular to the length axis. The 
upsetting operation is performed to increase the cross-sectional area of the material 
at the expense of its length. To achieve the upsetting, force is applied in a direction 
parallel to the length axis. 

There are four types of forging methods which are generally used. They are 

a. Smith forging - This is the traditional forging operation done openly or in open 
dies by a black smith, by manual hammering or by power hammers. 

b. Drop forging - This is the operation done in closed impression dies by means 
of the drop hammers. 
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c. Press forging - Similar to drop forging, the press forging is also done in dosed 
impression dies with the exception that the force is a continuous squeezing type applied 
by the hydraulic presses. 

d. Machine forging - Unlike the drop or press forging where the material is drawn 
out, in machine forging, the material is only upset to get the desired shape. 

During the application of the load in a forging operation, friction occurs at -the 
interface of tool and workpiece. This frictional loss causes an increase in the required 
load. This frictional loss also depends on the length of the workpiece. In open die up- 
setting, due to the friction at tool-workpiece interface, the free surface of the workpiece 
bulges or barrels and fracture may occur at barreled surface. 


1.1 Literature Survey 


The problem of forging of metals has theoretical as well as technological significance 
in that solutions would help to bring about a better understanding of the behaviour of 
metals during plastic deformation at high pressure. Such solutions would give specific 
answers to problems actually encountered in practice. It is for these reasons that 
a number of investigators have dealt with the forging problem and have proposed 
certain theoretical solutions which are based on a material having the properties of 
an ideal plastic solid. This leads to so-called slip line solution. Such solutions have 
been proposed by Prandtl [l], Hill, Lee and Tupper [2], and Green [3]. The slip line 
theory has been well developed to analyse nonhomogeneous plane strain deformation 
of a rigid-perfectly plastic isotropic solid. This theory can be used to give very good 
first-approximations to loads required to perform operations and provide indications 
of the manner in which material deforms. However, the work hardening effects can not 
be incorporated in this method and therefore a detailed solution for the deformation 
that takes place during continued loading beyond the yield point can not be obtained 
using this technique. Shabaik [4] investigated the bulging in upsetting by using the 
slip line theory. 

Hoffman and Sach [5] proposed the slab method. A complete analysis of the slab 
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method has been presented by Altan [6] for the axisymmetric closed die forging. In this 
method, a slab of infinitesimal thickness is selected and the equation representing the 
force equilibrium for the slab is derived, assuming that the deformation is homogeneous 
within the slab and that the major directions are the principal stress directions. The 
resulting differential equation is then solved with appropriate boundary conditions. 
The whole procedure is elaborate though the principle is simple . While this method 
can give very good predictions of the load variation with deformation, it is inherently 
incapable of predicting shape changes such as barreling in open die; forging. 

The upper bound theorem was formulated by Prager and Hodge [7]. If surfaces 
of velocity discontinuities are included [8], then it states that among all kinematically 
admissible velocity fields v*, the actual one minimizes the following expression : 

J* = -^= J v a eq ^0.5 e*j dV + j s t |An*| ds - ^ % v* ds 

Here e*j is the strain rate field derived from v* and | Av*| denotes the discontinuity in 
the tangential components of the velocity across the surface Sr- Further, a eq represents 
the equivalent stress, r denotes the shear stress on the surface of discontinuity 5p and 
Ti is the specified traction on the surface S t - The first term expresses the power spent 
in causing the internal deformation over the volume V of deforming body. The second 
term represents the power loss over the surfaces of velocity discontinuities including the 
boundary between tool and material. The last term represents the power supplied by 
specified tractions. What the upper bound theorem states is that the actual externally 
supplied power is never higher than that computed by using above equation. 

Kudo [9] applied the upper bound theorem to the problem of axisymmetric cold 
forging and extrusion. Here, he examined two types of admissible velocity fields. He 
introduced the conception of a unit cylindrical deforming region to simplify the analysis 
of complicated problems. Upper bound solutions for problems with Coulomb’s friction, 
non-steady processes and work hardening material are also included in this paper. Kudo 
[10] also applied the upper bound theorem to the analysis of the plane strain forging 
and extrusion. Avitzur [11], S. Kobavashi [12], McDonald, Kobayashi and Thomson 
[13], Kobayashi and Thomsen [14,15], Liu [16] and many others have suggested upper 
bound velocity fields to predict the forming load in upsetting with bulging. 
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The major difficulty in using the upper bound method is how to choose the 
kinematically admissible velocity field since the accuracy of the solution depends on 
how close the assumed velocity field is to the actual one. 

Accurate determination of forging parameters under realistic conditions became 
possible when the finite element method (FEM) was introduced. The major advantage 
of the finite element method is that the method can be applied to a wide class of 
boundary value problems without restrictions of workpiece geometry. This is achieved 
by the proper discretization procedure used in the finite element method. Using FEM, 
it is possible to predict the platen forces, the pressure distributions, the strain (i. e. 
metal flow ) and hardness distributions within the billet, for various deformation levels 
of interface shear stress. 

Early applications of the FEM to forging problems were based on the incremental 
method proposed by Lee and Kobayashi [17]. The method uses the elastic-plastic stress- 
strain matrix based on the Generalised Hooke’s law and the Prandtl-Reuss equations. 
The additivity of incremental elastic and plastic strains is assumed. Even though the 
stress-strain matrix and the geometry are updated after every increment, only the 
linearised incremental equations are used. The method was applied to solid cylinder" 
upsetting [18], ring compression [19] and for predicting defects in upsetting [20]. Lee 
and Kobayashi [21] studied the problem of simple upsetting with complete sticking at 
the tool-workpiece interface using the method. Maccarini et al [22] used the method for 
studying the influence of die geometry on cold extrusion forging operation. Whereas 
Lee and Kobayashi [17] used the velocity as the primary unknown; Hartley et al [23] 
proposed an incremental method with the displacement as the primary unknown. In 
their method also, the linearised incremental equations are used and the elastic-plastic 
matrix and the geometry are updated after every increment. They studied the the upset 
forging process using the method where the interfacial friction has been incorporated 
using the /3-stiffness method [24]. 

Linearised incremental equations give only an approximate solutions. If the incre- 
mental size is not sufficiently small, the error between the exact and the approximate 
solutions grows rapidly with the applied load, as the error in the solution of the cur- 
rent increment gets propogated into the next increment when the stress-strain matrix 
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and the geometry are updated. To avoid this phenomenon, the non-linear incremental 
equations should be solved by using a scheme like Newton-Raphson scheme. Such a 
formulation was first proposed by Bathe [25]. The updated Lagrangian formulation of 
Bathe [25] has been applied by Carter and Lee [26] to axisymmetric upsetting. The 
thrust of the work of Carter and Lee [26], however, has been towards studying the in- 
terface mechanics rather than analysing the axisymmetric upsetting. They have used 
the contact surface algorithm to enforce the friction condition (Coulomb’s law) at the 
interface. 

Even though a lot of finite element studies have been carried out on axisymmetric 
upsetting, most of them have used linearised incremental equations without using suffi- 
ciently small incremental size. As a result, these analyses are approximate . There is a 
need to carry out a more accurate analysis of the process using the updated Lagrangian 
formulation with a Newton-Raphson iterative scheme for solving the incremental equa- 
tions. This is the objective of the present study. Of course, for an accurate analysis, 
the interfacial friction also needs to be modelled properly. This requires development 
of an appropriate contact algorithm. However, due to lack of time, this part of the 
work could not be carried out. 


1.2 Modelling of Cold Forging Process 


In the present work, the process of cold upsetting of cylindrical solid disk is studied. 
The process is considered to be axisymmetric. The updated Lagrangian formulation 
which is convenient for handling geometric and material nonlinearities is used. Newton- 
Raphson iterative scheme is used for solving the incremental equations. The material is 
assumed to be elastic-plastic strain hardening yielding according to von Mises criterion. 
The effects of temperature and strain rate ( viscoplasticity effects ) on the yield strength 
of the material are ignored in this work. The inclusion of these effects renders the 
analysis quite complex. Due to small acceleration, inertial forces are not included. 
The body forces are also neglected. 

The dies are assumed to be rigid and rough. So no relative displacement is 
permitted at the tool-workpiece interface and the process is considered under condition 
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of complete sticking at the tool-workpiece interface. 

1.3 Objective 


In the present work, a finite element code for a large deformation elastic -plastic prol>- 
lem is developed using the updated Lagrangian formulation of Bathe [25]. The code is 
developed for displacement controlled problems. 

First, the code is validated by comparing the results with experimental as well 
as FEM results of reference [21]. Then, a detailed parametric study is carried out. to 
show the effects of three process variables namely the percentage reduction in height 
(% r), the height to diameter ratio (H/D) and the material properties. The param- 
eters considered are the forging load, the contact normal stress distribution at the 
tool-workpiece interface, the effective stress and strain contours and some deformation 
patterns. 

The deformation pattern is studied using the deformed mesh or the bulging pro- 
files of free surfaces of the cylindrical solid work- piece. The prediction of the metal 
flow in the workpiece is done by using the contours of the effective strain. 


1.4 Plan of the Thesis 

The thesis is organised as follows : 

In the second chapter, the updated Lagrangian method for problem with geomet- 
ric and material nonlinearities and the corresponding finite element formulation are 
presented. The application of the boundary conditions and the solution procedure for 
the problem of cold upsetting are discussed at the end. 

In chapter 3, the results obtained by applying the method for a range of cases 
are presented along with the comparison with experimental results and a discussion on 
the parametric study. 
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Conclusions and suggestions for further work are discussed in chapter 4. 



Chapter 2 


Mathematical Modelling of an 
Axisymmetric Upsetting Process 

In this chapter, the mathematical model for cold upsetting process is developed The 
process is modelled as a steady state axi-symmetric problem. The constitutive equation 
for the mechanical behavior of the work material is stated. The interaction between the 
work and it’s surrounding are represented by appropriate boundary conditions. The 
finite element method formulation is presented at the end of this chapter. 


2.1 Governing Equations 


In the study of the deformation of a body subjected to external loading, often the 
original undeformed and unstressed state of the body is used for the formulation of its 
equation of motion: This is known as Lagrangian Formulation. This formulation is con- 
venient for small deformation which is the case in a majority of engineering problems. 
In such cases, the deformed configuration does not deviate much from the original one 
and hence the deformation can be described by an infinitesimal strain tensor, known 
as the engineering strain tensor, for which the strain-displacement relation is linear. 

On the other hand, for large deformation problems one has to use a finite strain 
measure which is expressed by a non-linear strain-displacement relation. Furthermore, 
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the equations of motion of the reference configuration depend on the deformation. 
Hence, for all large deformation problems the Lagrangian Formulation proves to be 
cumbersome with the governing equations being difficult to solve. In such cases, one 
solves the problem using an incremental method known as Updated. Lagrangian Formu- 
lation. In this formulation, it is assumed that the state of the stress and deformation of 
the body is known till the current configuration, say at time t. The main objective is to 
determine the incremental deformation and stresses during an infinitesimal time stop 
d t, i. e. from time t to t + dt. Here, the current configuration is used as the reference 
configuration for obtaining the incremental values. Unlike the Lagrangian Formula- 
tion, an incremental strain tensor is used. The methodology is particularly useful 
for elastic-plastic materials because the stress-strain relationship in such materials is 
usually expressed in an incremental fashion. 

2.1.1 Strain-Displacement Relations 

Figure ( 2.1) shows a body in its current configuration (which is being used as the 
reference configuration) at time t and a possible deformed configuration at time t + At. 
The vectors and u t+At denote the displacements of a point P of the body at time 
t and i + At, respectively. The position vectors of P at the above times are depicted 
as ri and r t+At , respectively. For a very small value of At, one can represent it as dt. 
Thus , from the figure it can be seen that du is the incremental displacement of P for 
the time step dt. For small incremental deformation, one can represent it tensorially 
using the following linear incremental strain tensor: 

dey = — (du tJ - + diij'i) (2-1) 

Note that de^ is not the incremental engineering strain tensor but the incremental 
logarithmic or natural strain tensor. 

If the incremental deformation is large, the incremental strain-displacement rela- 
tionship has to be non-linear, i. e. : 
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Configuration 
( time = 0) 


Figure 2.1: Configuration 


The tensor d e v - is known as Green Lagrange strain tensor. Henceforth, it is assumed 
that the incremental displacement is small and thus the linear strain-displacement 
relation given by eq. (2.1) holds true. 

2.1.2 Constitutive Equations 

The plastic part of the incremental stress-strain relationship is obtained from the plastic 
potential of the material. For a material, yielding according to the von Mises criterion 
in an isotropic fashion, the plastic potential is given by [27], 


F(o, 3 ,p) = crater,,) - CT»(p) 


(2.3) 


Note that 


F = 0 


(2.4) 
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represents the yield criterion. The plastic potential F depends on the Cauchy stress 
tensor a l3 through it’s second invariant a eq called as equivalent stress and defined by 


' eq 


f 3 / t 
= ( 3 ^% 


(2-b) 


where a is the deviatoric part of a ty Further, F depends on the variable yield stress 
of the material, a y , through a hardening variable p. For the case of strain hardening, 
p is identified as the equivalent plastic strain e p eq , and hence defined as : 


P = e e«j = (2.G) 

and 

d«;, = (\<K,O k (2-7) 

Here, def., is the plastic part of the incremental strain tensor de tJ and the integration 
in eq.(2.6) is to be carried along the particle path.The dependence of a v on p ( or 
normally approximated by a power-law type of relationship 


= K(e p eq ) n ( 2 - 8 ) 

Here, K is called the hardening coefficient and n is called as the hardening exponent. 


The incremental plastic strain (de^) is obtained from the plastic potential using 
the following relation : 

FtF 

(2.9) 




u 


where dA is a scalar. This equation is called as the flow rule. Differentiation of eq.(2.3) 
with respect to < 7 ^ gives 


dF 


dc r. 


o 


2a ° tJ 


Using this and from eq.(2.9), one can determine dA as: 


dA = def. 


( 2 . 10 ) 


( 2 . 11 ) 


Further, the hardening relationship and the yield condition (eq. 2.4) can be used to 
express dA as : 

j \ _ dcr y _ da eq 


( 2 . 12 ) 
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where 


H __ 

d^’ 


( 2 . 13 ) 


the slope of the hardening curve, can be obtained from eq. (2.8). Substitution of eqs. 
(2.10) and (2.12) in eq. (2.9) leads to the following constitutive equation 


p 3 d(7 e g i 

d( " = 2 


( 2 . 14 ) 


This constitutive relationship between the deviatoric stress tensor and the in- 
cremental plastic strain tensor is not really convenient for an updated Lagrangian 
formulation for which the incremental stress-strain relationship is needed. This can be 
obtained from eq. (2.14) as follows: 


j p 3 O" d(T e q 

d<, = o' £r b - -do-fci 

1 2 H<j eq da kl 


( 2 . 15 ) 


Note that, from eqs. (2.3) and (2.10), we get 

doeq _ dF_ _ _3_ 
daki dcr kl 2cr e 


(2.10) 


Substitution of eq. (2.16) in eq. (2.15) leads to the following incremental plastic stress- 
strain relationship : 

d < - < 2 - 17 > 


The incremental elastic stress-strain relationship is given by 
de* 3 = i[-i/d cr kk 5ij + (1 + i/)d a tj ] 


( 2 . 18 ) 


where defj is the elastic part of de^ , E is the Young’s modulus and u is the Poisson’s 
ratio. Adding the two relationships, we get 


d eij = de-j + de? 


- [-&*** + t ** + iiw,\ Aau 


^ a ij a kl 


( 2 . 19 ) 
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This is the incremental elastic-plastic stress-strain relationship needed in the updated 
Lagrangian formulation. However, it is the following inverse relationship [28] which is 
more useful : 

da, 3 = 

where 

9 /^ [jor'ki 


C5Z = 2 /, 


5ik$ji + 


1-2 u 5ljSkl 2(3 n + H)al 


eqj 


and m is the shear modulus related to E and v by the relation: 

E 


M = 


2(1 + Z/) 


( 2 - 20 ) 

( 2 . 21 ) 


(2.22) 


In eq (2.20), if the increment consists of only a pure rigid rotation, then da t] must 
be zero. Such a stress increment is called an objective stress increment. Thus, in the 
incremental elastic-plastic stress-strain relationship (eq. 2.20), the incremental stress 
tensor must be an objective stress tensor. The Cauchy stress tensor does not satisfy t.ho 
objective requirement. So doy, in eq. (2.20) can not be interpreted as the incremental 
Cauchy stress-tensor. The two most commonly used objective stress tensors are (i) the 
2nd Piola-Kirchojf stress tensor (ii) stress tensor whose material time derivative 
is given by the Jaumann stress rate cr tJ . If we use the latter, then day,- in eq. (2.20) 
should be identified as oy, dt. Thus 

<jij dt = C^JdcH • (2.23) 

The Jaumann rate a \j is related to Cauchy rate a\j by the relation 

oy, dt = a tJ dt - cr lfc (0, fc dt) - a jk (Q lk dt) (2.24) 

where the incremental rotation tensor (Q tJ df) is defined as 

Q tj dt = ^(du M - d u Jft ). (2.25) 

Equations (2.21, 2.23-2.25) represents the elastic-plastic incremental stress-strain rela- 
tionship to be used in the updated Lagrangian formulation. 
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2.1.3 Equilibrium Equations 

In the upsetting, we assume that the load is applied slowly so that the resulting acceler- 
ations are very small. Then, the inertial forces become negligible and the problem can 
be treated as a static problem. Neglecting body forces also, the equilibrium equation 
becomes 

= 0- (2.26) 

This equation is to be satisfied in the deformed configuration (i. e. in the configuration 
at time t+dt) which is not known yet. In updated Lagrangian formulation, we need 
equilibrium equation in terms of the incremental stress. The incremental form of eq. 
(2.26) is given by [29]: 

dcqj ) j (<jjj j i)d'iifcj — 0. (2.2 1 ) 

Note that dcr tJ in the above equation is the incremental Cauchy stress tensor and hence 
equal to a tJ dt. 

Now the above equations can be collected together to define the updated La- 
grangian formulation. 

2.1.4 Updated Lagrangian Formulation 

It is assumed that the net load is applied in several increments and after a certain 
number of increments, we arrive at the current configuration (i. e. the configuration 
at the present time t), which is treated as the reference configuration. Further, the 
displacement vector u t and the stress tensor cr,j at the present time are known. Now 
a load increment is applied to the current configuration. Then to find the resulting 
incremental displacements, strains and stresses (du„ de^, doy,-), we have to solve the 
following equations: 

(a) Incremental equilibrium equations: 


dOijj ((Jtj^ldufcj — 0; 


(2.28) 
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(b) Incremental stress-strain relations: 

o-y d t = Cf : [ L de kl ; (2.29) 

[ oy, defined by eqs. (2.24-2.25) and C% kl defined by eq. (2.21) ] 

(c) Incremental strain-displacement relations: 

dfijj = —{du t j + duj tl ) (2.30) 

Note : When the incremental deformation is large, the non-linear incremental strain 
displacement relation (eq. 2.2) should be used. 


2.2 Boundary Conditions 


The typical boundary conditions are as follows: 

(i ) Geometric or essential boundary conditions. 

The incremental displacement vector is specified at a point on the boundary. 

(ii) Force or natural boundary conditions. 

The incremental stress vector d U — day, nj is specified at a point on the boundary, 
where n 3 represents the components of the unit outward normal vector. 

Sometimes, we have mixed boundary conditions, i. e. at a point on the bound- 
ary, some incremental displacement components and the remaining components of the 
incremental stress vector are specified. 

Figure 2.2 shows the domain of the problem. Due to the symmetry, only the 
upper right quarter of the work piece is used for analysis. 

The boundary conditions for the domain of our problem are discussed below: 

(i) The workpiece-die interface (AB): 

Because of assumption of complete sticking condition at the workpiece- die in- 
terface, the component of incremental displacement in r-direction is zero. Also, since 



2.2 Boundary Conditions 


16 


Z 




Uppei 

A 

platen 

B 


D 


C 


Lower platen 


Figure 2.2: Compression of a cylindrical solid disk between two platens 


our problem is displacement controlled, the incremental displacement component in 
the ^-direction is specified. 

Thus, the boundary conditions on AB are 

Essential boundary conditions: 

du r = 0, d u x = specified (2.31) 

(ii) The free surface ( BC ): 

Since this boundary is a free surface, it is traction free. When the workpiece 
comes in contact with the die , this surface bulges. Therefore, the components of 
incremental displacement du r and du 2 can not be assumed to be zero. Thus, the 
boundary conditions on BC are : 

Natural boundary conditions: 

d t z = dcr rz = 0 (2.32) 

and 

dt r = d CTrr = 0 
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(iii) The plane of symmetry (CD ): 

Because of symmetry, the z-component of incremental displacement and r-component 
of incremental stress vector are zero. So we have the following boundary conditions on 
CD: 


Natural Boundary Condition: 


d t r = d o rz = 0 

Essential Boundary Condition: 


(2.33) 


d u z = 0 (2.3-1) 

(iv) The axis of symmetry, z-axis (DA) : Again, because of symmetry, the 
r-component of incremental displacement and z-component of incremental stress vector 
are zero on DA. Therefore, the boundary conditions on DA are 

Natural Boundary Condition: 


d t z = d a rz = 0 


(2.35) 


Essential Boundary Condition: 


du r = 0 


(2.3G) 


2.3 Finite Element Formulation 


In this section, a finite element formulation corresponding to the updated Lagrangian 
approach of section 2.1 has been developed. First, an elastic-plastic matrix relating 
incremental stress and strain vectors is derived using the plastic potential proposed in 
section 2. 1. Then, the equilibrium equations in the deformed configuration are derived 
using the principle of virtual work following a procedure given in [25]. Finally, the 
numerical scheme required to implement the finite element formulation, is discussed. 
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2.3.1 Incremental Stress-Strain Relationship Representing 
Material Non-linearity 


Since the stress-strain relationship is path dependent for an elastic-plastic material, it 
is convenient to represent it in incremental form. The incremental stress-strain relation 
is expressed in vector form as 



{da} = [C\ EP {de} 

(2.37) 

where 


{da} T = {da rr ,dcr zz ,da rz ,dagg} 

(2.38) 

and 


{de} — {de rr , 2dc rz , 

(2.30) 


are respectively the vector forms of the incremental (objective) stress and strain tensors 
The expression for the elastic-plastic matrix [ C EP ] is derived from the plastic potential 
F. For its derivation, it is convenient to express the yield condition eqs. (2.3) and (2. 1) 
as 


F{{a},p) = cr e<7 ({cr}) - cr y {p) = 0 

(2.40) 

where 

{a} = {a rr , cr zz , cr rz , CJgg} 

(2.41) 

On the yield surface, dF = 0 and hence 


, 9F . r . , dF 

W (dcr} + = 0 

(2.42) 

or 

{a} T {da} — AdA = 0 

(2.43) 

where, the flow vector {a} is defined by 


r , dF 

° 

(2.44) 

the parameter A is given by 

„ 1 dF A 

A = dA9p dp ' 

(2.45) 
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and dA is the same scalar which appears in eq. (2.9). Eq. (2.10) gives the following 
expression for the flow vector: 

W = J-M'. (2.4C) 

where {cr}’ is the vector form of the deviatoric part of oy,-. Eqs. (2.6) and (2.11) give 


dp = dA = de p eq . 


(2 47) 


Substitution of eqs. (2.40) and (2.47) in eq. (2.45) leads to the following expression for 


A: 


A = d°y _ dcr v 
dp de? 9 


( 2 . 48 ) 


The total strain increment can be split into elastic and plastic parts. Thus 

OF 


{de} = {de 6 } + {de p } = [Cp^a} + dA 


d{a) 


= [C] ^dcrj + dAla} 


( 2 . 49 ) 


Here, [C] is the matrix form of the elasticity tensor C tJ ki (See eq. 2.18). The second 
part of the right hand of eq. (2.49) is due to the flow rule (eq. 2.9). Premultiplying 
both sides of equation (2.49) by {o} T [C] and using eq. (2.43) , we get 


> M r P]W 

AL + {oHC]{a} 


(2.50) 


Substituting the above expression for dA in eq. (2.49), we get 



Ml rn-WH-rl 1 M T [ C K<i£} 

{de} — [C] (M + A + { a yr[ C ]{ a }M 

(2.51) 

which leads to . _ . . . . . 


r j i [r] [CIWW'P] mi 

{dcr} - \[C] A + ^ a |T[ C j { a } J 

(2.52) 

Therefore, 


rri *P ( lC \ \ 

[C] "V 1 A + {a}T[C\{a}) 

(2.53) 


where the expression for {a} is given by eq. (2.46). For a material with power strain 
hardening law (eq. 2.8), the expression (2.48) for A reduces to 


A = Kn(e p q r- 1 


(2.54) 
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For an isotropic material, the expression for [C] 

for the 

axisyinmotiu 

• case is given by 


' i - V 

V 

0 

V 


fci — e 

(1 + i/)(l 2v) 

V 

0 

1 -u 

0 

0 

\-7v 

2 

V 

0 

(2.55) 


V 

V 

0 1 

-i/ j 


As stated at the end of subsection 1.1.2, the vector {dcr} 

in eq, 

. (2.52) should be 


interpreted as an objective stress increment. 

2.3.2 Updated Lagrangian Formulation for Geometric Non- 
Linearity 

When a body experiences a large deformation and/or rotation, the equilibrium must, 
be established in the current configuration. The normal solution procedure adopted in 
this case is updated Lagrangian. That is to say that when all the kinematical variables 
are known from time 0 to time t in discreate time steps, the objective is to establish 
the equilibrium in the configuration at time t+At. In order to derive the finite element 
equilibrium equations at time t + At , the principle of virtual work requires that 

j ^ cr£ At S(e^ At ) dV t+At = R t+AL ( 2 . 56 ) 

where cr^ At is the Cauchy stress tensor at time t + At , 8(e\^ AL ) is the virtual infinites- 
imal strain tensor at time t + At , R t+At is the virtual work of the surface forces at 
time t + At (body forces being neglected) and V t+At is the domain at time t + At . 
The major difficulty in application of eq. (2.56) is that the configuration at t + At is 
unknown. Morever quantities like Cauchy stress tensor are not purely additive as one 
has to take care of rotation into account. 

An elegant way of formulating the problem is given by Bathe [25]. Here, the 
virtual work expression (eq. 2.56) is transformed to an integral over the known domain 
at time t (i. e. V t ) using the 2nd Piola- Kirchoff stress tensor and the Green-Lagrange 
strain tensor. Thus, the virtual work expression becomes 

/ ,S‘* Al S(,e‘t A ‘) dV* = # +i ‘. 


(2.57) 
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where 5 tJ denotes the 2nd Piola- Kirchoff stress tensor whose relation with the Cauchy 
stress tensor a ZJ is given by 

* 5 i +ai = (2.58) 

[The right superscript stands for the current configuration and the left subscript denotes 
the reference configuration.] 


Since the Cauchy stress tensor is always referred to the current configuration, the 
left subscript has been omitted. Here f+At^m represents the derivative of the position 
vector r l at time t with respect to the one ( r t+At ) at time t + A t and -^rst denotes 
the ratio of densities at time t and 1 4- At. The Green-Lagrarige strain tensor ,e[* Ai is 
defined as 


, ,.C-KA£ I ,.,C+iA£ t+L 

* e t] — 2 [ iU tJ +tU j,* + t u *, i tu kd 


(2-59) 


where t u j+ At is the derivative of the displacement vector u t+Ai at time t + At with 
respect to the position vector at time t. For axisymmetric case, the strain component 
in 0-direction is given as 


A ie Af = 


U t+At + i / n t+At ' 


An incremental decomposition of stress and strain gives 


where, 

^for axisymmetric case 
and 


tS^ Al = tSl, + tAS lJ 

= + tASxj , 

5 t e£ At = S t Ae XJ + S t At?,.,, 


tAcij — — (tAtitj “b tAuj'i) 

,A £M = ^ ) 


tAriij — — ( t Auk,i) (tAufcj) 


(2.60) 

(2.61) 

(2.62) 


(2.63) 



2.3 Finite Element Formulation 


22 


(for axisymmetric case 




Here t Au stands for the derivative of t A u ( incremental displacement vector at time 
t ) with respect to r 4 (position vector at time t ). Thus, eq. (2.57) with incremental 
decomposition can be written as 


J vt tAS l] 5( t Ae t j)dV t + t A5 I /( t A7 7lJ )dy t + J^a^Av^dV 1 

+ J aljSit A€ t3 )dV l = 7? 4+At . (2.64) 

The 2nd integral on the left hand side is a higher order term compared to other terms 
and hence can be neglected. Moreover, using an elastic-plastic incremental stress- 
strain relationship (eq. 2.20), jASV, can be approximated as L C t Ae k i- Here, the 
left subscript of denotes the time at which it is to be evaluated. The above 

simplification will definitely give rise to an error in the right hand side of eq. (2.64). 
The error is normally neutralised by using some iterative scheme like Newton Raphson 
Method which is discussed in subsection 2.3.4. Now, eq. (2.64) can be written as 

J yi t C,% t Ae u S( t Ae,,)dV‘ + A^dV + J^SUA^W 

= (2.65) 


Owing to the symmetries in tC^ii tAc y , tAr) XJ and cr^, eq. (2.65) can be rewritten 
in the following form: 

f vt (i,{Ac} T ) t {C\ EP ,{Ae}dV‘ + J vi (5[{Atj} t ) [T]‘,{A>j}<1V‘ 


m 


t _ 


+ J vi (<S,{Ac} T ) {<r}'dV‘ = 

(2.6G) 

-rr> t Ac zz? 2 t Ae rz , iAcqq} , 

(2.67) 

{a^r, o\ z , <x‘ z , <7go} T 

> 

(2.68) 

t&Ur^z, t&Uz^r, t&U Zi z<, J ’ 

(2.69) 

El* 0 0 ■ 

o E1‘ 0 
o o El* . 


(2.70) 


where, 
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and 



* CTrr 

at T Z 

0 ‘ 

H‘ = 

°ZT 

<4 

0 

0 

0 

<4 . 


(2 71) 


The t [C] matrix is the elastic-plastic matrix evaluated at time t. which is given by 
eq. (2.53). 


2.3.3 FEM Equations 

The domain is discretised into a number of elements and the incremental displacement 
field over a typical element e is approximated as 

({M = [Q] t{Au} e (2.72) 


where, 

£ {Au} e = { f Au>, t Auj, . . . , t Au", ( A<}. 

and t Au* , t Au* stand for the (unknown) incremental displacements of node Y in r 
and z directions respectively . Here 


{Qi} T 

m T 

where 



(2.73) 


{Qi} T = [ Ml 0 JV, 0 ... AT. o], 

{«W r -[0 JV, 0 N 2 ... 0 JV„]. (2.74) 

and N t , the functions of ( r, z ), are called the shape functions. 


With displacement field defined by eq. (2.72), the strain field within the element 
can be calculated in terms of nodal displacements. From eq. (2.72), 

• A ^ = ^ =,Wi} 5‘ {A " } ‘ (2 ' 5) 


m \T __ 


3{Qi}’ 


to 


where 
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Substituting eq. (2.75) in eqs. (2.62) and (2.63) and arranging them in vector form, we 
get the strain displacement relations: 

t {Ae} = t [B L ] t{A«} e , (2.77) 


t {A77} = t [B N ] t {Au}' 

where, 

tmi 
t{Q*% 

tm T r + t{Ql} T z 

As il r 

T t 


,mi 


and 

.[Bivf = [ .{(30J 


(2.7S) 


(2.70) 


(2.80) 


Using the strain displacement relations derived above, the contribution to the 
integral (eq. 2.66) from a typical element e (with volume V* ) can be written as 

MA«}' r t[ cr ,{Bl] dV‘ j ,{ Au}' + 

5 t {A K }' T ^ t [B M f [T)‘ ,[S„ ] ,{Au} e + (2.81) 

MAu}‘ T ^ t [B L f {<r}'dvA = <5,{A U }‘ r ( t+ a.W) 

Here, the contribution to the term R t+ from the element e has been expressed in terms 
of the elemental external force vector t+At{F} e using a standard procedure. Since the 
variation in the incremental displacement vector 5 t {<5u} e is arbitrary, expressing the 
term within 1st parenthesis as linear elemental stiffness [ Ki ] e and that of 2nd one as 
nonlinear elemental stiffness [ K^l ] e the above equation can be written as 

( t [K L ] e + t [K NL f ) t{Au} e + t {f} e = t+A *{F} e (2.82) 


or 




(2.S3) 
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where 

*{/}' 6 = J yt t[B L ] T {<r} £ dVA (2.84) 

is the internal force vector. 

Assembling the elemental matrices t[A] e and the elemental vectors t{/} e and 
(+At{^} e over all the elements, we get the following global equation : 

t[K] ({Au} + t {/} = t+At{F}. (2.85) 

Here t[F] is the global stiffness matrix, ;{Au} is the global incremental displacement 
vector (at time t ) and t{f} and t+&t{F} are the global internal and external force 
vectors respectively. Decomposing t+At {F}, eq. (2.85) can be written as 

t[K] t{Au} + t {f} — t {F} + t {A.F}. (2.80) 

Here, t {F} is the (global) external force vector at time t and t {A F} is the (global) 
incremental force vector from time t to t + At. In updated Lagrangian formulation, it 
is assumed that the equilibrium equations have been satisfied exactly at time t. Thus, 

t{/} = t{F}. (2.87) 

Then, eq. (2.86) reduces to 

t [K\ t {A u} = t {AF}. (2.88) 

Now, suppose we find the incremental stresses corresponding to the solution 
( t {Au}) of eq. (2.88), add those to the earlier stresses to find the stresses at time 
t + At and use eq. (2.84) to determine the internal force vector t+ A i{/} at time t + At, 
it is possible that it may not match with the applied (or external ) force vector. Thus, 
eq. (2.88) represents only an approximate equilibrium equation at time t+ At (The ap- 
proximation is mostly due to the linearization and simplification involved in the steps 
between eqs. 2.64 to 2.65). A solution of such an approximate equation may involve a 
significant amount of error and, depending on the time/load step, may becomes unsta- 
ble. Therefore, it is necessary to modify eq. (2.88) to turn it into an iterative problem 
capable of providing a solution with desirable accuracy. The details of the iteration 
scheme is discussed in next section. 
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2.3.4 Numerical Scheme 

The finite element equation of an elastic-plastic material (eq. 2.88) can be solved it- 
eratively by various numerical schemes. Newton Raphson scheme is one of the most 
commonly used method for such problems. However, if the incremental load step is too 
large, the error involved in the computation of stress from the [C]^^ 5 matrix, evaluated 
at the previous point, also becomes large [30]. In such a situation, the backward return 
algorithm can be used to minimize the error. Another limitation of Newton Raphson 
method is that it often fails to converge in the neighborhood of critical point. 

Normally we encounter two types of problems in mechanics: (1) Load control 
problem and (2) Displacement control problem. In a load control problem, the desired 
deformation is achieved by prescribing a load at a point. On the other hand, in a 
displacement control problem, it is the prescribed displacement which gives the desired 
deformation. In a load control problem, if the load falls due to change in geometry, then 
most iterative schemes fail to converge. Displacement control problems are free from 
this limitation. In this thesis, only displacement control problems are considered. 

(I) The Newton Raphson Scheme 

The scheme is best described with the help of Fig. (2.3) in the following steps. In 
Fig. (2.3), m is an equilibrium point at load level {P} . The aim is to reach another 
equilibrium point ( m + 1 ) at a load level ( {P} + A{P} ). 

Step 1 — First of all, the stiffness matrix t [K] in eq. (2.88) is evaluated using the 
stresses corresponding to point m. Next, eq. (2.88) is solved with *{AP} = {AP} t.o 
get an incremental displacement vector t {Aw}. 

Step 2 — From the incremental nodal displacements, incremental strain compo- 
nents are calculated by using eq. (2.77). Incremental stresses are calculated from 
the incremental strains by using suitable constitutive relation. Generally incremental 
stresses are calculated at the Gauss points of the elements. If the point under con- 
sideration is in elastic range, the elastic stress-strain matrix [C] is used. Otherwise, 
elastic-plastic matrix [C\ EP is used. Since elastic-plastic matrix is dependent on the 
current state of stress, the stress value at position m is used for the computation of 
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Figure 2.3: Illustration of Newton Raphson iterative scheme for a single load 

elastic-plastic matrix. Eq. (2.23) is used to calculate Jaumann stress rate and the total 
Cauchy stress follows from eq. (2.24). 

Step 3 — Once the stresses are calculated at all the Gauss points of all the ele- 
ments, the internal nodal force vector t+At{/} is calculated. This force corresponds to 
the equilibrium load at point i in Fig. (2.3). 

Step 4 — Since the applied nodal forces ( t {P} + t {AP} are not equal to internal 
nodal forces t+At {/}, the deformed configuration is not an equilibrium configuration at 
load level (*{P}+ t {AP}). It gives rise to an unbalanced nodal forces ( t {P} + *{AP} - 

t+At{/})- 
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S_ tep 5 — The steps 1 to 3 are repeated with unbalanced nodal forces as the load 
vector and i as a starting point. Successive repetition of the above steps will lead 
to the equilibrium point (m + 1) for load level t {P} -f t {AP} when the unbalanced 
nodal forces become almost zero(less than a preassigned quantity). In each successive 
iteration incremental displacement vector is added to its previous value. 

Faster convergence is achieved if the tangential stiffness matrix at each successive 
equilibrium points is updated. The general practise is to achieve the total load in several 
increments. Within each incremental load step, Newton Raphson iterative scheme is 
applied with the same tangential stiffness matrix. Stiffness matrix is updated only 
when going for next load increment. This modified Newton Raphson scheme is a good 
compromise between reduced computational time and faster convergence. 

(II) The Backward Return Algorithm 

Consider the increment from point m to point i (See fig. 2.4). In calculation of the 
stresses at point i ( in step 2 of the Newton Raphson scheme ), we use the incremental 
stress-strain relation (eq. 2.37) where the [C] EP matrix is evaluated at point m. When 
the incremental load step is too large, this gives stresses not at point i but at some 
point i which is not on the yield surface (see Fig. 2.4). Thus 

M, = Mm + [C]* P { Ae}. (2.89) 

A further return is required to bring the stresses to point i. Therefore we can write 

Mi = M,' + {M- ( 2 - 9 °) 


To determine {Acr}, we proceed as follows. From eq. (2.50) and (2.53) we note 

that 

[C\ E m P {Ae} = [C] {Ae} - dA m [C]{a} m . (2.91) 

Therefore, to bring the stresses to point i, a correction is needed only in the second 
term. We assume that the correction can be approximated as 

{A<r} = -dA„[C]{a}... (2.92) 



Strain 

Figure 2.4: Illustration of Backward return algorithm 


To determine the scalar dA 0 , we expand the plastic potential at point x . Thus, 

* = + {Aff} + l^ Ap+ --' (2 ' 93) 

Analogous to eq. (2.45), we define A { > as 

, 1 dF t , . , 

Substituting eq.s (2.92) and (2.94) in (2.93) and using the definition of flow vector {a},/ 
(eq. 2.44), we get 

Fi = F t > - dA c {a}T [C]{a} t » - dA 0 A { > . (2.95) 

Since the plastic potential has a zero value on the yield surface, F x = 0. Thus 

dA 

° {a}J[C\{a}< + A< 


( 2 . 90 ) 
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Substitution of eq.s (2.92) and (2.96) in (2.90) gives the following approximate expul- 
sion for stress at point i: 




FACKah 


{a}'([C]{a} i / + A t 


A further expansion will lead to the result 


(2.97) 


{a}i » {<r} t » 


K'lCM' 


FiACKa}'. 


{«.}^[C]{a} t ' + A t > {“}/» [C] {«}*" + A t " 


(2.98) 


where the point i is closer to point i than i (see Fig. 2.4). When the plastic potential 
at point i becomes less than a pre assigned small number, no additional terms need to 
be added to equation (2.98). 


2.4 Choice of Shape Functions and Numerical In- 
tegration 


Since the finite element method is a numerical procedure for solving complex engineer- 
ing problems, important considerations pertain to the convergence of the numerical 
solution. For monotonic convergence, the elements must be complete and the elements 
and, mesh must be compatible. If these conditions are fulfilled, the accuracy of the 
solution increases continuously as we continue to refine the finite element mesh. This 
mesh refinement should be performed by subdividing a previously used element info 
two or more elements; thus, the old mesh must be ” embedded ” in the new mesh. 
The requirement of completeness of an element means that the polynomial approx- 
imating the displacements of the element must be able to represent the rigid body 
displacements and the constant strain fields. The requirement of compatibility means 
that the displacement within the element and across the element boundaries must be 
continuous. 

Eq. (2.65) suggests that the convergence criterion will be satisfied if du r and du* 
are chosen to be bilinear functions of r and z. To satisfy the above condition, we 
have chosen the second order Serendipity approximation with eight noded rectangular 
elements for the finite element formulation. Here, the term r' 2 z 2 of the assumed ap- 
proximation polynomial is omitted because it increases the size anti the bandwidth of 
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the stiffness matrix without improving the appi oximation outside the element. The 
approximation is normally expressed in terms of the natural coordinates £ and 77 . Thus, 
the shape functions in eq. (2.74) are the second order Serendipity shape functions 
in £ and rj. 

Exact evaluation of integrals appearing in element coefficient matrices and source 
vectors is not always possible because of the algebraic complexity of the coefficients in 
differential equations. In such cases, it is natural to seek numerical evaluation of these 
integral expression. Numerical evaluation of integrals, called numerical integration, 
involves approximation of the integrand by a polynomial of sufficient degree, because 
the integral of a polynomial can be evaluated exactly. 

The most used numerical integiation method is Gauss-Legendre numerical in- 
tegration and formulae for integrals defined over a rectangular master element iln is 
given by 

J Qr *!)<%&! = J ^ F(£>v) d^j d£ (2.99) 

m n 

~ £ £ F {£uVj)uJiWj ( 2 . 100 ) 

1=1 t=l 

where m and n denote the number of quadrature points in the £ and i] directions, (£„ r/j) 
denote the Gauss points, and w t and iuj denote the corresponding Gauss weights. 


2.5 Evaluation of the Secondary Variables 

As explained above, solution of the problem is obtained in the form of nodal displace- 
ments. From nodal displacements, the secondary quantities like the upsetting load, the 
normal stress distribution at the tool-workpiece interface, the contours equivalent stress 
and strain and the deformed mesh are calculated. The details of these calculations are 
explained briefly in the following lines- 

i. Upsetting Load : 

Since complete sticking is assumed at the die- workpiece interface, there is no 
frictional loss at the interface. Then the supplied energy is exactly equal to the work 
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spent in creating the internal deformation. Thus, over the increment from t to t + At 
( t AF z )( t Au z ) = cr'+ A< t A( eq 27rr f+At dr t+At dz' tAt (2.101) 

where tAu z = specified platen displacement in the increment 


a, 


t AF z = required incremental upsetting load 

/3 tt+At /t+At\ 2 

T = ) . 

= equivalent stress based on stress tensor at time t + 6t 

/2 \ 3 

tA€ e q — l~tACij (ACfj J , 

= incremental equivalent strain from t to t + At 

t 

The incremental upsetting load is found from eq. (2.101). 

ii. Surface Contact Stress Distribution : 


( 2 . 102 ) 


(2.103) 


The surface contact stress distribution is evaluated by considering each of the 
billet elements in contact with the platen and by linearly extrapolating from the values 
of a zz at the Guass points. 

iii. Equivalent Stress and Strain : 

The equivalent stress ( a eq ) and the equivalent strain ( e eq ) are calculated at each 
Guass point from the components of stress and incremental strain tensors at that Guass 
point by using eq. s(2.102) and (2.103) and the following relation : 

t+At 

e et A = E 

Then, the contours of a cq and e eq are plotted using the values evaluated at 2x2 Guass 
points of the elements. 

iv. Deformed Mesh : 

The deformed mesh of the workpiece is plotted to exhibit the deformation pattern. 
The deformed mesh after any increment can be obtained by updating the r and z 
coordinates of nodal points by using the displacement components upto that increment. 



Chapter 3 


Results And Discussion 


The finite element modelling of an axisymmetric upsetting process developed in the 
previous chapter has been applied to a number of cases involving two metals and var- 
ious sets of input variables to illustrate its applicability. The metals considered are 
aluminium HOOF and steel AISI 1019 whose material properties are given as below. 
The properties for the aluminium are taken from [21] and that of the steel from [23]. 

Material properties : 

Aluminium 1100F 

Young’s modulus E = 69:11 GPa 
Poisson’s ratio v = 0:33 
Intial yield stress Y 0 = 89:85 MPa 

Slope of effective stress and plastic strain curve h' = 138:23 MPa 
Steel AISI 1019 

Young’s modulus E = 210 GPa 
Poisson’s ratio u = 0:3 
Intial yield stress Y 0 = 478:2 MPa 
Hardening parameter K = 563:98 MPa 
Hardening exponent n = 0:14 
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In this chapter, first a convergence study is performed to check the convergence of 
the finite element formulation with respect to the mesh size. Then, the finite element 
code is validated by comparing the results with the FEM and experimental results 
obtained by Lee and Kobayashi [21]. In the last section, the parametric study is carried 
out to illustrate the effects of various process variables on the process parameters. 


3.1 Convergence Study 

The convergence study is performed for the problem of an axisymmetric upsetting of 
an aluminium 1100 F cylinder with the height to diameter ratio of 1. 2. As stated in 
chapter 2, only a quadrant is chosen for the analysis. Initially, the mesh of 35-eight 
noded rectangular elements is selected (fig. 3.1). The elements along the die contact 
interface, along the free surface and along the axis of symmetry are made smaller as 



the stress/displacement gradients are likely to be high in these regions. Then, the 
refinement is performed by dividing each element into four elements forming a mesh of 
140 elements. The results are obtained for both the meshes. 
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Figure 3.2: Convergence of load-displacement, curve. 


Fig. (3.2) shows the variation of non-dimensional load ( the ratio of the upsetting 
load (L) to the initial yield strength of the metal ) with percent reduction in height for 
the two cases. It can be seen from the figure that as we refine the mesh, the solution 
seems to converge. There is a decrease of 8.89 % in the load at 20 % reduction when 
the mesh is changed from the coarse mesh to the fine. It is also found during this study 
that a smaller increment of die displacement resulted in a slightly lower value of the 
upsetting load. Thus, the solution seems to converge as the sizes of both the elements 
and the displacement increment are reduced. 

Fig. (3.3) shows the convergence of the bulge profile, i. e. the displacement of the 
free surface of the cylinder at 20 % reduction in height for the two meshes. 

Figures (3.4a) and (3.4b) show the contours of equivalent strain at 20 % reduction 
in height for 35 elements and 140 elements respectively. It is found that the values of 
equivalent strain also have converged. 
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r/Ro 

Figure 3.3: Convergence of bulge profile at 20 percent reduction in height. 


INCREMENTAL STEP = 0.03 
MINIMUM VALUE = 0.03 



INCREMENTAL STEP = 0.02985 
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MAXIMUM VALUE = 0.5373 



Figure 3.4: Convergence of contours of effective strain at 20% reduction in height 
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3.2 Validation 

The computer code is validated by comparing the results with the FEM and experimen- 
tal results of Lee and Kobayashi [21], They performed the analysis of an axisymmetric 
upsetting process with the height to diameter ratio of 0.8 using a mesh of 124 four 
noded quadrilateral elements for half the cylindrical specimen. The material proper- 
ties of A1 1100F were used for the analysis which was carried out upto a 33 percent 
reduction in height. The full sticking condition was used at the die-workpiece interface. 

The same problem is solved using the developed finite element code. A quarter 
of the cylinder specimen is analysed using the mesh of 140 eight noded rectangular 
elements with smaller size elements along the die contact area, along the free surface 
and along the axis of symmetry. 

Fig. (3.5) shows the comparison of the FEM and experimental load-displacement 
curves. The discrepancy between the theoretical and experimental predictions can 



Figure 3.5: Comparison of theoretical and experimental load-displacement curves. 

be explained as follows. In experimental studies, slip is observed over a part of the 
interface [21]. Whereas in theoretical modelling, full sticking is assumed over the whole 
interface. As a result, the theoretically predicted behaviour of the elements on the free 
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surface near the edge shows that they undergo a very large deformation ( see figures 
3.6a and 3.7a ), which doesn’t happen in real life upsetting process. Consequently, 
the theoretical predicted values of contact stress ( near the edge of free surface ) 
become very' high as shown in fig. (3.9). Even though the theoretical model neglect the 
interfacial slip, the energy needed to create the above mentioned large deformation of 
the free surface elements seem to be more than what is dissipated due to the interfacial 
slip. This is the reason why the theoretically predicted values of the upsetting load are 
higher than the experimental values. The values of Lee and Kobayashi [21] are lower 
than the values predicted by the present study. This is because the deformation gets 
under-predicted ( see fig. 3.7b) by the method of Lee and Kobayashi as they used the 
linearised incremental equations. 


This phenomenon of very large deformation of the free surface corner elements is 
observed when H/D < 1 or %r > 20. Therefore, these cases are avoided in most of the 
parametric study carried out in section 3.3. 





Figure 3.6: Bulge profiles at various reductions in height. 


Figures (3.6a) and (3.6b) show the comparison of the bulge profiles at various 
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reductions in height. For small reductions the theory well predicts the geometrical 
change of the free surface. The predictions, however, deviate from the observation as 
the reduction in height increases. After some stage of compression, the phenomenon 
of free-surface fold-over is predicted at the corners. The extent of this relative upward 
metal flow at the corners increases with friction and with the amount of deformation. 
In real upsetting, flow past the dies is not possible. Therefore, the originally free surface 
comes into contact with the die at the corners [31]. 



INCREMENTAL STEP = 0.03133 
MINIMUM VALUE = 0.03133 
MAXIMUM VALUE = 0.46995 



Figure 3.7: Distribution of effective strain at 20 percent reduction in height (a) present 
study (b) results of [21]. 
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The effective strain is an indication of the degree of deformation, and can be 
calculated by following the deformation at any point incremently. The comparison of 
the effective strain distributions at a 20 percent reduction is shown in fig. (3.7a) and 
3-'b). A percentage error of 11.2 is observed with respect to the values of Lee and 
Kobayashi [21]. 


3.3 Parametric Study 


In this section, a study is made of the effect of process variables viz. the percent 
reduction in height ( %r ), the height to diameter ratio ( H/D ) and the material 
properties on the important design parameters like the upsetting load, the contact 
pressure distribution at the interface, the bulge profile and the distributions of effective 
strain and effective stress. 

3.3.1 Effects of Percent Reduction in Height and Height to 
Diameter Ratio 

In this section, analysis is carried out for three cases of height to diameter ratio of 
cylindrical specimen ( H/D = 1, 1.25 and 1.5 ) and for various values of percentage 
reduction from 5 to 20 percents. The material properties of aluminium HOOF are used. 

1. The Upsetting Load 

Fig. (3.8 ) shows the relation between the upsetting load and the percentage 
reduction in height (i. e. the die displacement). The relation is linear in the 
elastic range as well as in the early stage of yielding, but it deviates after the 
plastic region spreads considerably. The upsetting load increases with reduction 
in height. Besides, as the figure indicates it decreases with H/D ratio. In real 
materials in which stress-strain curves are not straight but always convex to the 
upper direction, the load decreases with H/D ratio [32]. 
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Figure 3.8: Load-displacement curves for various H/D ratios. 


2. Distribution of Contact Pressure on the Interface 

Figures (3.9a and 3.9b) depict the distributions of normal contact stress on the 
interface of die and workpiece at 5 and 20 percent reductions in height respectively 
at three different values of H/D namely 1, 1.25, 1.5. The normal stress values 
increase with radius, r and tend to infinity at the edge r = R 0 . This happens 
because of the assumption of full sticking condition. Similar trend is reported by 
Lee and Kobayashi [21]. If the interfacial slip is taken into account, the normal 
stress value is likely to be finite at the edge. 

The normal stress distribution doesn’t change either with H/D or with percent 
reduction. However, the values increase with percent reduction ( for the same 
H/D ratio ) and decrease with H/D ratio (for the same percent reduction ). 

3. Bulge Profile 

Figures (3.10a and 3.10b ) show the bulge profiles at 5 and 20 percent reductions 
respectively at three values of H/D (H/D=l, 1.25 and 1.5). Figures indicate 
that as upsetting advances, the bulge becomes larger. It is because the elements 
near the edge, which yield early suffer a large deformation with the increase in 
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(a) 



(b) 

Figure 3.9: Distribution of normal contact stress at interface for various H/D ratios at 
( a ) 5 percent reduction ( b ) 20 percent reduction. 


reduction and the growth of bulge depends mainly on the deformation of these 
elements. At the same percentage reduction in height, the size of bulge increases 
in proportion to the H/D ratio. 
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4. Distribution of Effective Strain 

Figures (3.11a), (3.11b) and (3.11c) show the contours of effective strain for the 
following three cases respectively : (a) %r = 20, H/D = 1, (b) %r = 20, H/D = 
1.5 and (c) %r = 33, H/D = 1. There is no significant change in the distribution 
pattern of effective strain when either %r or H/D is changed. However, the values 


INCREMENTAL STEP = 0.02932 INCREMENTAL STEP = 0.03023 

MINIMUM VALUE = 0.02932 MINIMUM VALUE = 0.03023 



f I fLo 

(c) 


Figure 3.11: Contours of effective strain ( a ) %r = 20, H/D - 1 ( b ) %r - 20, H/D 
= 1.5 and ( c ) %r = 33, H/D = 1. 

of effective strain decrease with H/D ratio and increase with %r. 



3.3 Parametric Study 


45 


5. Distribution of Effective Stress 

Figures (3.12a), (3.12b) and (3.12c) show the contours of effective stress for the 
following three cases respectively : (a) %r =20, H/D = 1, (b) %r = 20, H/D = 
1.5 and (c) %r = 33, H/D = 1. 



In the early stages of upsetting, the regions at the center and at the edge of the die- 
workpiece interface yield first ( > 89:85 MPa). As upsetting advances, these 
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regions spread and join each other. The figures show that, at 20% reduction, the 
region at the center of the interface has still not yielded. A small region near the 
free surface also has low values of equivalent stress as this region is un-restrained. 
At 33% reduction, the whole region becomes plastic. 

The figures show that as H/D is increased, the values of effective stress decrease 
and the distribution pattern also changes. 

6. Deformed Mesh 

Figures (3.13a), (3.13b) and (3.13c) show the deformed mesh for the following 
cases : (a) %r = 20, H/D = 1, (b) %r = 20, H/D = 1.5 and (c) %r = 33, H/D 



Figure 3.13: Deformed mesh ( a ) %r - 20, H/D - 1 ( b ) %r 20, H/D 
( c ) %r = 33, H/D = 1. 
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1. The specimen with H/D ratio of 1.5 undergoes more deformation than 
that of H/D 1. Further, it is observed that upto 20% reduction in height, the 
phenomenon of free-surface fold-over doesn’t appear as long as H/D ratio is larger 
than 1 . However, whenever %r increases beyond 20 or H/D becomes less than 1 . 
the fold-over of the free surface occurs as shown in fig. (3.13c). 

3.3.2 Effect of Strain Hardening 

In this section, the process parameters are obtained for two different materials namely 
aluminium 1100F and steel AISI 1019 to study the effect of strain hardening on them. 

1. The Upsetting Load 



Figure 3.14: Load-displacement curves for aluminium and steel (H/D — 1.25) 

Fig. (3.14) gives the comparison of load-displacement curves for steel and alu- 
minium. The upsetting load required to deform the steel specimen is obviously 
larger than the one required for the aluminium specimen as the yield stress is 
larger. However, the load-displacement curves are similar even though the hard- 
ening laws are different. 
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3.3 Parametric Study 


2. Distribution of Normal Contact Stress at the Interface 

Distribution of normal contact stress at tool-workpiece interface for the two met- 
als is shown in fig. (3.15). The normal stress for steel is larger than the one 
for aluminium over the whole contact surface. This again is due to the larger 
values of the yield stress of steel. There is not much change in the contact stress 
distribution pattern of the two metals. 



Figure 3.15: Distribution of normal contact stress at interface for two metals ( %r = 

20, H/D = 1.25) 


3. Bulge Profile 

Fig. (3.16) shows the bulge profiles of steel and aluminium specimens at 20 per- 
cent reduction in height and for H/D = 1.25. It is observed that bulge profiles 
are more or less the same. This shows that the bulge profile depends purely on 
the kinematic parameters like %r and H/D ratio and not an strain-hardening. 

4. Distribution of Effective Strain and Effective stress 

Fig. (3.17 ) shows the distributions of effective strain for aluminium and steel 
while fig. (3.18) shows the contours of effective stress for the two metals. The 
figure show that the pattern and the values of the effective strain are almost 
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the same for the two metals. This implies that the effective strain distribution, 
like bulge profile, depends only on the kinematic parameters. The distribution 
patterns of the effective stress are also similar. The values for the steel specimen 
are of course larger. 
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Figure 3.17: Contours of effective strain at 20%r and H/D — 1.25 for ( a ) aluminiu 
and ( b ) steel. 



Chapter 4 


Summary, Conclusions and Scope 
for Further Work 


In the present work, the axisymmetric upsetting problem under the condition of full 
sticking at the die-workpiece interface is solved by the finite element method considering 
the material to be elastic-plastic strain hardening yielding according to von Mises 
criterion. The inertial terms are neglected. The updated Lagrangian formulation which 
is convenient for handling material and geometric nonlinearities is used. The modified 
Newton- Raphson iterative scheme is used for solving the incremental equations. Since 
the method permits to follow the deformation in step-by-step manner, it is possible to 
determine the quantities at any instant. The usefulness of the method for carrying out 
detailed studies of the deformation and stress distributions in plasticity problems has 
been demonstrated. 

The detailed parametric study has been carried out choosing aluminium as a 
material to manifest the effects of the process variables like the percent reduction in 
height (%r) and the height to diameter ratio (H/D). The important observations are as 
follows. It is obvious that the upsetting load and the values of the normal contact stress 
increase with %r. However, the'pattern of the normal contact stress distribution along 
the interface doesn’t change much with %r. The pattern does not change much with a 
change in H/D ratio also. However, the values of normal contact stress decrease with 
H/D ratio and, as a result, the upsetting load also decreases but only to some extent. 
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Along the interface, the normal stress increase with the radial coordinate. When a 
full sticking condition is used, the normal stress tends to infinity at the edge. This, of 
course, is not realistic. Therefore for better modelling of the process, the interfacial 
slip should be accounted for. 

The size of the plastic zone is determined using the contours of effective stress. 
The plastic zone first appears at the center and at the edge of the die-workpiece in- 
terface. With increase in %r, these two regions spread and join. But even at as high 
value of %r as 20, the region at the center of the interface remains elastic. Only at very 
high %r, the full region becomes plastic. The size of the plastic zone reduces with H/D 
ratio but only by a small amount. As for as the pattern of effective stress distribution 
is concerned, it changes with H/D ratio but not with %r. 

The deformation pattern is studied with the help of deformed mesh and the 
contours of the effective strain. The bulge of the free surface becomes larger when 
either %r or H/D ratio is increased. However, the effect of H/D is small. When a 
full sticking condition is used, the phenomenon of the fold-over of the free surface is 
observed for %r > 20 or H/D < 1. The pattern of effective strain distribution doesn’t 
change much either with %r or H/D ratio. However, the values of effective strain 
increase with %r and decrease with H/D ratio. 

The influence of strain hardening is analysed by obtaining the results for two 
materials namely aluminium and steel. The values of the normal contact stress, the 
upsetting load and the effective stress, of course, larger for steel than for aluminium. 
However, the patterns of the normal contact stress distribution along the interface, 
the variation of upsetting load with % r and the effective stress distribution are not 
different. It implies that strain hardening has no effect on these three parameters. It is 
observed that the bulge profile as well as the values and distribution pattern of effective 
strain also do not depend on strain hardening. 

The important conclusions are : 

• Along the interface, the normal stress increases with radial coordinate. When 

a full sticking condition is used, the normal stress at the edge acquires a non- 
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realistic value of infinity. Therefore, for realistic modelling of the process, inter- 
facial slip should be incorporated. 

• When a full sticking condition is used, the phenomenon of the fold-over of the 
free surface is observed for %r > 20 or H/D < 1. This, again, doesn’t happen in 
a real-life upsetting process. 

• The H/D ratio doesn’t seem to have much influence on the process parameters. 

• The deformation (i. e. the bulge profile and the values and the distribution pattern 
of effective strain doesn’t get affected by strain hardening. 

• Strain hardening has no effect on the pattern of normal contact stress distribution, 
the variation of upsetting load with %r and the effective stress distribution. 

The following works may be taken up as extensions of the present study : 

1. Proper modelling of interfacial friction and slip. 

2. Dependence of strain rate (viscoplasticity) and temperature on material behavior. 

3. Analysis and prediction of defects in the upsetting process. 

4. Warm upsetting process. 
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